Sub-Analysis: Variation of Model Constants for a Circumbinary System¶

A 488 Computational Astrophysics Project by Tyler Barna¹

¹Department of Physics and Astronomy, Rutgers, the State University of New Jersey, New Brunswick

Dependencies¶

For my sub-analysis, I use almost all of the same libraries as the main analysis. However, I import some libraries slightly differently and provide different aliases that I find convenient. Additionally, my $\texttt{plotter.py}$ is slightly altered from the main notebook; I added labels for the specific bodies as well as units for the axis. The only additional library I import is $\texttt{warnings}$ Which I use to surpress an expected error later in my analysis.

3-Body Problem¶

Since I'm only planning on using a 3 body simulation and I'd like to vary two different constants, I'm going to create a function that takes the parameters for Kepler-693 (The subject of Jonathan's Sub-Analysis) in the case that the stars form a Circumbinary Orbit, which makes it possible to estimate the initial positions of all three bodies using Kepler's Laws. This is accomplished by using a reduced mass for the first star and planet to determine the posittion of the second star. In effect, this method partially simplifies the system to a psuedo-two body system in order to determine the iniital conditions for all three bodies. Since I'm not changing any of these, the only arguments I have for this function are the number of orbital periods to graph, the value of the Universal Gravitational Constant, the smoothing length, and which graph I want to be used from $\texttt{plotter.py}$

The evaluation in this cell has the same parameters that Jonathon uses in his analysis and is what I will be using as my "control" to compare against when I vary the constants.

Varying the Universal Gravitational Constant¶

Despite traditionally being accepted as a constant value, various studies since 1962 have derived values of G that have disparities greater than their uncertainties. According to a 2015 paper by Anderson et. al, it is possible that G varies with time periodcially, with an estimated amplitude of $1.619 \pm 0.103 \times 10^{-14} \: \mathrm{ m^{3} kg^{-1} s^{-2} }$. This has potentially novel implications for the stability of this circumbinary orbit.

Minimum Value of G¶

Maximum Value of G¶

Varying the Universal Gravitational Constant: Conclusions¶

Looking at each of the evaluations, it seems that there isn't a particularly striking difference between these and the control system. This isn't wholly surprising, as the range of values for G is only around ~$0.02$. It's also worth noting that the validity of this study has been called into question, with a 2015 commentary by M. Pitkin suggesting that there is most likely an additional noise term that explains the perceived periodicity.

Varying the Smoothing Length¶

The smoothing length is perhaps a more salient constant to vary, as it is a constant chosen largely at our own discretion for this project and introduces a degree of bias to the simulation. I have outlined several different cases where I vary the value of the smoothing length from the original value of $\epsilon_{orig} = 0.001$

Case 1: $\epsilon = 0$ ($\epsilon_{orig}\times 0$)¶

Setting $\epsilon = 0$ breaks the model, which makes sense; $\epsilon$ is used to prevent numierical divergences. I have supressed the error warnings here since they are quite long and were anticipated.

Case 2: $\epsilon = 0.00001$ ($\frac{ \epsilon_{orig} }{100}$)¶

Case 3: $\epsilon = 0.0001$ ($\frac{ \epsilon_{orig} }{10}$)¶

Case 4: $\epsilon = 0.01$ ($\epsilon_{orig}\times 10$)¶

Case 5: $\epsilon = 0.05$ ($\epsilon_{orig}\times 50$)¶

Case 6: $\epsilon = 0.1$ ($\epsilon_{orig}\times 100$)¶

Case 7: $\epsilon = 1.0$ ($\epsilon_{orig}\times 1000$)¶

Varying the Smoothing Length: Conclusions¶

There proved to be quite a significant difference in simulations with smoothing lengths equal to a multiple of the original value. While each of the cases has a fairly unique outcome, I will attempt to draw some hollistic conclusions from the overall results

For smaller values of $\epsilon$:¶

As the smoothing length gets closer to zero, the system increasingly begins to behave similarly to a two-body system, as the stars begin to have a more standardized orbit with one another and behave similarly to a single star. Very low values are not entirely realistic, as in a physical system the stars would not be simple point masses and would have a certain minimum distance from one another to prevent collision.

For larger values of $\epsilon$:¶

As the smoothing length increases, the systems become more and more erratic, with a value somewhere between $\epsilon=0.05$ and $\epsilon=0.1$ being the point that circumbinary orbit is no longer stable for this system, with $\epsilon=0.01$ being possibly unstable. While it is important to at least evaluate these extreme values, there is a point where the smoothing length is so large that it no longer accurately portrays a realistic system and the simulation stops being informative.

Overall, I wasn't able to determine a particularly consistent way of determining an appropriate smoothing constant value besides through trial-and-error. However, I suspect there are theorists much more experienced than myself who are able to determine less-biased heuristic methods using statistical analysis.